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I. Introduction. 

Spectral element methods are high order methods that combine the flexibility of finite 
element methods with the "infinite order accuracy" of spectral methods. The domain of computation 
is decomposed into some subdomains - the elements - (generally these are deformed 
parallelotopes) , and the exact solution is approximated by a piecewise polynomial of high degree. 
The spectral element method differs from other spectral methods using domain decomposition 
techniques (patching methods) by the way the matching conditions are handled. These are, like in 
the finite element method, implicitly taken into account by the variationnal statement of the 
discrete problem. This allows for more flexibility with no loss of the spectral accuracy ( see, e.g. , 
[P], [M.P.] and [F]). When the algebraic equations resulting from this kind of discretization are 
obtained, the problem that remains is to solve, in an efficient way, the algebraic system. 

The interest of domain decomposition technique is to fraction the computational task so as to 
yield smaller problems and to use parallel computers for instance. If the value of the approximate 
solution were known on the various interfaces, the problem would be very simple since it would 
results into the resolution of as many disconnected problems as the number of elements. The main 
difficulty is that this value is not known; hence a technique known as the iteration per subdomains 
has been proposed in the literature [F.Q.Z.] [Q.S1] to discover this value iteratively, Another 
approach is to try to invert the whole problem by not working iteratively on each subdomain. This 
method, [M.P.] [F.R.D.P.] [R], consists in reducing iteratively the residue at the same time on 
every subdomain. The global method used can be based on a conjugate gradient algorithm or another 
iterative procedure. 

In the first part of this paper, E. M. Ronquist and A. T. Patera [R.P.] have presented some 
results concerning a new multigrid method for the resolution of the algebraic system resulting 
from the approximation of a second order P.D.E. by spectral or spectral element method in the 
one-dimensional domains. The very simple idea of using the Jacobi preconditioner as a smoother 
for the iterative multigrid algorithm appears to be a very good one. Indeed, the numerical 
properties of this smoother seems to surpass all expectations; the reduction rate of each V-cycle 
appears numerically to be independent of the discretization parameters. 
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When iterative techniques are used, it is important to understand why these methods 
converge in order to foresee the generalization and the ability of the methods to be adapted to more 
than test problems. Here we propose an analysis of this phenomenon and provide the justification 
of these very good properties. Many general convergence proofs exist in the literature for the 
numerical analysis of the multigrid technique; among them let us cite [M.M.], [B.D.], We use here 
the abstract framework developed by R. E. Bank and C. C. Douglas [B.D.] that fits exactly the 
numerical conclusion of [R.P.] concerning the optimal choice of the smoothing operations. 

To our knowledge, the numerical analysis of the convergence for the multigrid algorithms 
used in spectral type techniques is somehow empty. The main reason is certainly that the nice 
analysis that can be done requires a variational framework and the awareness that the spectral 
methods are, exactly or very close to, variational approach is not so old. The other reason, 
perhaps, is that the previous multigrid techniques applied to spectral type methods [Z.W.H. 1 ,2] 
used a finite difference preconditioner as a smoother and a Chebyshev framework. The 
convergence, in this case, is not so brillant as in the present approach and a priori more related to 
the good properties of preconditioner of the finite difference operator. Besides the variational 
formulation involves a nonsymmetric form that makes the analysis much more involved. 

The paper is organized as follows; in section II, we recall the theory of [B.D.] in a form 
adapted to our analysis. In section 1 1 1 , we first explain on the test example of the Galerkin spectral 
approximation of the homogeneous Poisson problem the fundamental reasons of optimal properties 
of this multigrid method. The tools are based on Jackson inequality and some refined version of the 
approximation properties of the L 2 -projection operator. In section IV , we generalize the analysis 
to the case of the spectral element approximation. We compare in each section the results obtained 
by the theory with the numerical results presented in [R.P.], The last section V deals with the one 
domain multigrid technique when applied to a non constant second order problem. 

The generalization of these results to multidimensional problems will be presented in a 


future paper. 
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II. Position of the problem and abstract theorem. 

11.1 Generalities on variational multigrid techniques. 

In this subsection, we first recall the theory developed by R, Bank and C. Douglas to analyze 
the convergence rate of the multigrid algorithm for solving the linear algebraic system that arises 
from the numerical approximation of elliptic partial differential equations. We present it in a 
version that we shall use afterwards. First of all, let 56 be a Hilbert space, a be a continuous, 
elliptic, symmetric, bilinear form and g be a continuous linear form, both defined over 56. The 
problem to be solved is : Find u € 56 such that 
(II. 1 ) V v € 56 , a(u,v) □ g(v) . 

For the numerical resolution of this problem , we first introduce a sequence of finite dimensional 
subspaces A, c cM 2 c ... c <M,- of 36 ; then consider the problem : Find u= € <M. such that 
( 1 1 ■ 2 ) j V v € cK j , a(Uj,v) = g(v) . 

The basic idea for solving problem (11.2) with a multigrid algorithm consists in first 
defining a simple problem over the largest spaces <M, . and solving it, then correcting the residual 
derived from the solution of this simpler problem when plugged into problem ( 1 1 . 2 ) j by solving 
problem (I l.2) k for lower values of k < j. The first step is the most important one and relies on 
the good choice of continuous, elliptic, symmetric, bilinear form b, called smoother, that 
represents a in some sense and is easier to invert. Let us suppose that we have only two grids, the 
coarse one (i.e., ) and the fine one ( i.e. , cH , 2 , and j = 2), The two-grid procedure consists in 

1 - m/2 steps of smoothing where we solve m/2 times a problem like the following one : Find 

c><p in oK 2 such that 

(11.3) V v € A 2 , b($(p-ip,v) ng(v) - a((p,v) , 

2- one step of coarse grid correction where we solve only once a problem like the following one: 
Find ip in such that 

(11.4) V v € , a(?,v) = g(v) - a(<p,v) , 

and define 0(p a <p + «p . 

3- m/2 steps of smoothing as in (11,3). 

We consider here just two grid levels, the only reason is for sake of simplicity of the 
notations, but as in [B.D.], we could consider the whole W-cycle$ based on more than two levels. If 
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the initial guess for the exact solution u is u° , after one V-cycle like the one described previously 
by the three steps 1 - 2- 3-, the resulting solution is u 1 and can be expressed like a function of u° 
as follows 

(n.5) u 1 = s m/2 es m/ V) , 

so that, after the r lh V-cycles, the solution is 

(11.6) u r = [S m/2 e£ m/2 ] r (u°) . 

Moreover, it is very simple to note that if u is the exact solution, then 
6u = u and <£u = u , 

so that, if e r denotes the error after the r lh V-cycle, we derive from (11.6) that 

(11.7) e r = [S m/2 eS m ' 2 ] r (e 0 ) . 

Note that the equations (11.3) (11.4) define affine operators 6 and § but in (11.7) we can consider 
these operators as linear ones since they operate on differences e° and e r ; from now on, we shall 
consider these operator as linear ones while keeping the same notation. 

As noted in [B.D.], of importance is the analysis of the spectrum of the following eigenvalue 
problem : Find ^ in cM 2 and X in R + such that 

(11.8) V v € cH 2 , a(4 ; ,v) = X b(^ ,v) , 

where b has been scaled so that the maximum eigenvalue X max is equal to 1. Let us order the 
eigenvalues in increasing order 0 < X, < X 2 < ... < X p = X max = 1 , (where P is the dimension of 
cH- 2 ) and choose relative eigenfunctions 4 ; 1 , ^ 2 .... .St'p . Of equal importance in the analysis is 
also the compatibility between the coarse space ^ and the space spanned by the first 
eigenvectors . 

More precisely, under the following hypothesis 

H the space l, coincides with span {4' 1 , tP 2 , ...^ }, (where p is the dimension of cA l, ) ; 
then one can prove the following theorem 

Theorem 11.1 : Under the hypothesis H the error after the first V-cycle verifies 
a(e 1 ,e 1 )<(l-X p+1 ) 2m a(e°,e°). 


Proof : First of all , let us recall that the eigenfunctions , n = 1 .... ,P , form a basis of l 2 that 
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is orthogonal for both the forms a and b. Let us span e° in this basis; we get 
e° = L n N =1 e n ° ; 

then, due to (11.3) and (11.8) , we derive that c> m/2 (e°) satisfies 
(11.9) ^ m/2 (e 0 ) = L n N =1 (1-X n ) m/2 e> n ■ 

From hypothesis H, we then derive that the operator© truncates the previous spectrum so that 

e^ /2 (e°) = Z: n N =p+1 (1-X n ) m/2 e>„ i 

then the smoothing procedure diminishes once more the spectrum of the error as follows 
e 1 = $ m/2 0£ m/2 (e°) = L n N =p+ , ( 1 -X n ) m e n ° . 

We deduce now from the orthogonality of the that 


,e') = L " , ( I -X„) 2m ( e °) 2 at V „) < ( 1 -X p „) 2 ”L n N =P ., < 


\2m / ~0\2 


2m v* N 


^ 0\2 


< < 1 -Vi> 2 "T„ N s i ( S n 0 ) 2 a( Vn ) = ( 1 -X P+1 ) 2 " a(e°,e°) 


p+i 

p+i 


Remark 11.1 : Note that the previous theorem is very simple and is a trivial extension of the 
analysis of the multigrid procedure in the Fourier space. Note also that if the space cM 1 is not so 
well chosen, for instance if it coincides with span {SF N _ p+ } , '1 ; N _ p+ 2 , ..,^ N }, then the multigrid 
procedure would not converge rapidly since, after the first V-cycle, we obtain 
§ m/2 es rn/2 (e°) = L n N :? ( 1 -X n ) m e n ° , 

and the error remains important since X, can be very small. In fact, the method has exactly the 
same properties as the plain Jacobi algorithm. 


It turns out from the previous analysis that the multigrid procedure, when applied under 
hypothesis H, is degenerated since only one Y-cycle is needed and m is the only important factor of 
convergence. So in the non trivial applications where the hypothesis H is not verified, we have to 
measure the actual situation between the hypothesis H and the situation explained in remark 11.1. 
This can be explained as follow : the "rough" eigenmodes ( [B.D.] denotes this way the with X n 
close to 1) are damped during the smoothing procedure ( their components are multiplied by a 
factor ( 1 -X n ) ) while the "smooth" modes ( corresponding to small X n ) remain almost constant. 
Under the hypothesis H these are completely erased during the correction procedure, but if we are 
not in this optimal situation, they can only be damped also during this step. 
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The general analysis proposed in [B.D] allows for measuring the position with respect to 
both hypotheses. Let us recall the basis of their proof since we shall extend it in section V. To this 
purpose, they introduce the various norms, defined for any real index 0, as follows 

(11.10) III *111, = x, 8 } 2 ]’'' 2 . 

They introduce also the function f 

(11.1 1 ) f( a,b) = a a b b (a+b)' (a+b) = sup xe[0 ,j ( 1 -x) a x b . 

Then we derive that 

(11.12) Vim cM 2 , Vi , O<i<0 , HI £ m/2 <p ||| 9 - [Li, (1-X„) m X n %; ] 1/2 

< /(m/2,(0-l)/2) HI <p HI,. ; 

let us write now 

III es m/2 . ( in 2 □ a( e§ m/2 « , er ,2 . ( ) = a < es ,n,2 f . e§ m/2 * - r t > = .,< asr' 2 ., , r i,2 f) 

« ill es m,2 f in, _ 6 in s'" 22 .piii,., 

< [su P . e ^ 2 in e* in,.,/ in e* in, i in es" /2 » in, m s'" 22 * i« u , . 

so that we derive 

(n.t3) in es" 22 * iii, « [sup, „ M 2 in e* in,.,/ in e* m, j in s- 22 * m,., , 

which is valid for any £ > 0; we deduce from (11.12) that 

in s m/2 0s m/2 4, in, < f{ m/2, e/2) in e^ m/2 cp w M 

m/2, $/2) in es m/2 -p in, [sup^^jn 041 w 1 _/iii e* m, ] 

< K m/2.6/2) [sup, 6 ^ 2 III 041 HI,./ HI 041 1||, ] 2 ||| S m/2 4> ||| 1+p 

^ K m|) [sup, 6eAt2 HI 041 III,./ Ill 0^|||, ] 2 |||<p||| 1 . 

Defining^, 1 as follows 

A, 1 = { (p € cK 2 , such that V 41 eA, , a(<p ,41) = 0 } , 

( that coincides with the range of 0 ) and by minimizing the right-hand side over £, they state the 
following 

Theorem 11. 1 : Assume there exists a constant k ^ 1 and a > 0 such that for any <p € A, 1 

III * III,., < K* /2 III * III, . 

then 

(11.14) a(e' ,e’) < * a(e°,e°) , 
where 
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(11.15) -jf = [(k- 1 )/ K ] 2m if m<<x(K-1) 

* = (KV(m ,oc)) 2 if m>oc(K-1). 

11.2 Formulation of the spectral element discretization. 

Let us turn now to the position of the problem. We consider here the simple test problem 
over the interval A □ ]- 1 , 1 [ : Find u such that 

(11.16) - u xx = f , over A 

provided with homogeneous Dirichlet boundary conditions 

(11.17) u(- 1 ) = u( 1 ) □ 0. 

where f is a given force. This problem is very simple, but it allows a statement of the basic 
features of the multigrid algorithm and an understanding of why the method works. 

The spectral element method for approximating the solution of (11.16) consists in 
discretizing the space of acceptable functions by a subspace of piecewise polynomials. More 
precisely, given a pair h □ (K,N), we first break the interval A into K disjoint subintervals of 
comparable sizes 

A =U kr1 A k , A k =]a k ,a k+ b k [; 

then, we choose for space of approximate functions a subspace X^ 1 of Hg(A) consisting of all 
piecewise polynomials of degree < N, 

(II. 18a) X h N = Y h N fl H’(A) , 
where 

(11.1 8b) = { <p such that <p )Ak € IP N (A k ) } 

and P N (A k ) denotes the space of all polynomials of degree < N on A k . We remind that contrarily to 
the finite element approximation, the convergence is achieved by increasing the degree of the 
polynomials N and not refining the mesh. 

The discrete problem starts from the variational formulation of problem (11.16)01.17) 
that is : Find u in Hg(A) such that 
(11.19) V v € Hg(A), a(u,v) □ (f,v) , 

where (.,.) denotes the L 2 (A)-scalar product and a(.,.) denotes the following bilinear form defined 
over Hg(A) as follows 
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V ip,!|f 6 Hq(A) , a((p,«|r) □ (cp x ,^ x ) • 

Then we construct the numerical scheme by discretizing with Gauss-Lobatto quadrature formulas 
the various integrals present in (11.19) and restricting the space of test functions to X^ . 

This results in the following problem : Find u h € X^ , such that 

(11.20) V v h € \ , a h,GL^ U h > V h) “ ^’ V hVGL ’ 
where the discrete forms are defined as follows 

V (p ,4/ € Y h , (<p .t)h,GL = ^k=1 Z^ n= 0 ?n - 

V'P.'MXJJ , a h N , GL (‘P,»i<) ° (f x . >F x ) h % L • 

Here, the , and the are the weights and nodes of the Gauss Lobatto-Legendre formula with 
N+ 1 points and the collocation points L nk are defined by E; nk = a k + (E;^ + 1 )b k /2 . W recall here 

N 

that the integration formula is exact on P 2N _ 1 so that, contrary to the appearances, a h CL does not 

N 

depend on h nor N since a h GL = a ( see e.g. [M.P.]). 

The algebraic system that has to be solved is derived by choosing the values of the unknown 
function u h on the collocation points and representing u h in the basis of the interpolant basis h kn 
defined by : h kn is the only element of Y^ 1 such that 

<".20 = . 

The matrix system that has to be solved can be written as follows: Find u = (u£' n ) such that 

(11.22) A u = g , 

where A is the stiffness matrix with entries equal to 

(11.23) <2/b k )LT =l L" s0 cy k [dh„/dx dh jm /dx] U nl ) , 

I 

with p nk = P n b k /2 and denotes the direct "stiffness summation" , while g is related to the 
forcing term. We refer to [R.P . ] so as to [M.P.] for more details on the derivation of this matrix. 

This numerical method is proved to work very well in a great number of interesting 
problems that include, for instance, the full Navier-Stokes problem (see for instance [M.P.], 
[M.P.R.1], [M.P.R.2] or [R]) and is numericaly competitive and implementable on a parallel 
medium grain paradigm (see for instance [F.R.D.P.]). 


Let us turn now to the multigrid algorithm for solving iteratively problem (11.22). As 
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explained in [R.P.] , nested spaces are related to spaces of polynomials with lower degree X^ 1 (say 
M □ N/2) and the smoother is simply the Jacobi preconditioner B that is proportional to the 
diagonal part of the matrix A and normalized in such a way that the highest eigenvalue of B -1 A is 1 . 
We shall analyze in the two following sections the properties of this preconditioner and explain 
why the numerical method works so well as presented in [R.P.] . 

In this paper, we shall use two grids only though the analysis can be performed with no 
extra difficulty than the one of comprehension, and we shall use 
A, = X h N/2 and cK 2 a X h N . 

The problem enters in the general theory of [B.D.] since the numerical problem involves bilinear 
forms and matrices that do not depend on h. 


- 10 - 


III. Analusis of the convergence of the multiarid algorithm in the 
case of one .gJemenL 


As announced in the title, we shall assume for the moment that the discretization is applied 
on only one domain to problem (11,1 6). For sake of simplicity, we shall drop out the second index 
that corresponds to the subelement characterization in the notations. The situation is simple here 
as soon as we have explained the properties of the diagonal matrix, but this simple example 
permits enhancement of the main features that allow for a rapid convergence of the algorithm. The 
problem can be written as a pure collocation scheme as already noted in [M.P.], Indeed, by taking 
v h = h n for n = 1 ,N— 1 , in (11.20) and using twice the exactness of the Gauss-Lobatto formula, 
we derive 

a hVK’ h n) = a < U h’ h n) = - 1 a u h< x ) h nW dx = " ( U h , h n ) h ,GL = " U W 1 
besides we note that 

= e„. 

so that the problem actually verified by u h is : Find u h in X^ such that 
V n, n= 1 ,.. ,N- 1 , - u;'U n ) = fU n ) . 


This consists in a pure collocation procedure to solve the initial problem. 

In order to analyze the multigrid algorithm , let us first compute exactly the diagonal part of 
the stiffness matrix; that is here simply to the matrix with entries equal to 
\ = Ll_ 0 e n h;U n N )hjU n N ) , for 1,j = 1 N-1 , 

A Oj = A Nj n0, 


Due to the exactness of the Gauss Lobatto formula, we easily obtain that, for i = 1 N- 1 , 

A n = Ja [ h iM] 2 dx = - J A h i"M h i^ x ) dx ■ 

and here again, by using the exactness of the Gauss Lobatto formula, we derive from (11.21) 

As already derived in [GHO, formula (7,4)] 

Lemma III. 1 : The diagonal of the stiffness matrix verifies 
(111.1) A ii = Qj N (N+1)/[3(1-U, N ) 2 )] . 

Proof : Let us drop out in this proof, the superscript N. First of all, we note that, from the 
definition (11.21) of h. , we have 
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n":,' (x-y(x 2 -o 

h,<x) = _Jr: =«,n ; B :,'(x- Sj )(i-x 2 ) 

n":,' (E r V(5f-0 ^ 

j* i 

where a. is a non zero constant. It is well known (see, e.g [D.R.]) that the internal nodes of the 
Gauss Lobatto fomula verify 

(1 1 1.2) (x) = c N TTjlJ . 

where c N is a non zero constant. Hence, we can write 

(x-^j) hj(x) = c^ 1 (XjCl-x 2 ) L^(x) = 0Cj(1-x 2 ) L^(x), 
and after taking the derivative of both sides, we obtain 
(x-£.) h.’(x) + hj(x) □a i (d/dx)[(1-x 2 ) L^(x)] . 

Let us recall now the eigenfunctions property verified by the Legendre polynomials (see, e.g. 
[D.R.]) 

(111. 3) (d/dx)[( 1 -x 2 ) Lfl(x)] □ -(N)(N + 1 ) L N (x) ; 
we derive that 

(111. 4) (x-Kj) h.'(x) + hj(x) = -oc.N (N+1) L n (x) ; 
plugging now x = in this equation yields 

(111. 5) 1 n-oc,N(N + 1)L N (g . 

Besides, by taking the derivative of (1 1 1.4) , we obtain 

(1 1 1.6) (x-£.) hj"(x) + 2 h.'(x) o - a. N (N+1) L^(x) , 
plugging also here x = i;. in this equality and using (III. 2), we get 

(111. 7) hj'Uj) = 0 . 

Let us multiply now (1 1 1.6) by ( 1 -x 2 ) and take the derivative of the resulting equation, we deduce 
(x-Cj) (d/dx) [( 1 -x 2 ) h|"(x)] + 3( 1 -x 2 ) h."(x) - 4xh.'(x) = oc. N 2 (N+1 ) 2 L N (x) . 

Finally, plugging one more time x = , we derive from (III. 5) 

3 (1-£ 2 ) h.-U.) = - N (N+1) , 
and the lemma follows from the value of A.. . 

Remark 1 1 1. 1 : Note that, as a consequence of (1 1 1.7), we have proved here that the Lagrangian at 
point t- has its maximum at . 


Associated with the matrix A is the bilinear form a, and in the same way, associated with B , 
the normalized diagonal of A can be defined a bilinear form b ^ 1 CL .This will be the smoother of the 
multigrid algorithm. Lemma 1 1 1. 1 proves that the smoother we have introduced is proportional to 
the bilinear form b^ GL defined for any ip and »p in X^ (which here is simply 1P N ( A) fl H G (A)) 

(in. s) B" a (*,« l ,)=L" = „Piii(Ei)»<V(i-s?r'. 

It is interesting to note that, due to the exactness of the Gauss-Lobatto formula, we can rewrite 
b^ 0L in a continuous way since for any ip and 41 in X^ , the function [(p vp ( 1 -x 2 )~ 1 ] is still a 
polynomial and belongs to P 2 N _ 2 ^) . so that we have 

(111. 9) bJ GL ((p,ip) = B(<p,»p) = | A (p(x) xp(x) (1 -x 2 )- 1 dx . 

Let us now analyze the eigenvalue problem (11.8) or more precisely the eigenvalue problem 
associated to B. The situation is here very simple since we have an exact expression for the 
solution to this problem 

Lemma III. 2 : Let us define for any integer n , 1 < n < N- 1 

(111.10) ^ n (x) q(1-x 2 ) l;(x) ; 
then we have 

(111.11) V v € Xj , a(^ n ,v) = n(n+1 ) b(^ n ,v) . 

Proof : Let us first remind the following standard formula on the Legendre polynomial (see, e.g. 
[D.R., Chap 2, §7] ) 

(111.12) V n € IN , (( 1 -x 2 ) L;(x))' + n(n+l) L n (x) □ 0. 

Let us compute, for any v in X^ , 

a(^ n ,v) = J A v;(x) v'(x) dx = J A (( 1 -x 2 ) L;(x))' v’(x) dx ; 
using ( 111 . 12 ) and integrating by parts, we obtain 

a(4 / n ,v) = - n(n+ 1 ) j A L n (x) v'(x) dx = n(n+ 1 ) } A L' (x) v(x) dx 

= n(n+ 1 ) j A ( 1 -x 2 ) L„(x) v(x) ( 1 -x 2 )' 1 dx = n(n+ 1 ) b(^ n , v); 
this ends the proof of Lemma 1 1 1.2 . 

It is important to note that in this simple example the eigenvalues are well known and, 
moreover, that the first M of them span exactly the space IP M+ 1 (A) fl H G ( A). As a first 
consequence, we can state that the normalized form that will be used as a smoother is defined as 
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follows 

(111.13) Vip,H»€Xj , b((p,^) = N(N-1) b((p,q/) = N(N-1)J A tf(x) q^(x)/(1-x 2 ) dx . 

We demonstrate the case where the simple hypothesis H of section 1 1 is valid, and we can state now 
Theorem 111,1 : Let u° denote the initial guess in the multigrid algorithm applied to problem 
(11.20) in the case of one element, and u 1 denotes the solution obtained in the after m smoothing 
and one correction. This solution converges to the exact solution u as follows, 
aCu-u 1 , u-u 1 ) = [1 - (N + 2)/4(N- 1 )] 2m a(u-u°, u-u°) . 

P roof ; This is a simple corollary of Theorem 11.1 and lemma 1 1 1 .2 since the first ( N/2) — 1 
eigenvectors span the coarse space <Al 1 = X^ /2 ; it follows that the eigenvalues of problem 
(11.8) are n(n+ 1 )/(N- 1 )N, 

Remark II 1.2 : Note here that the correction on the coarse grid needs only be done once and that 
there is no optimal choice for the number of smoothings since the convergence is proportional to 
this number. This is actually in accordance with the numerical simulation of [R.P.] as appears in 
table 1 of that paper. 

Remark III. 3 : Let us point out the fundamental reasons that give this rapidity to the algorithm. 
They are hidden here due to the simplicity of the eigenvalue problem. First of all, even if this is 
not of major importance, the matrix B is a good preconditioner of the matrix A. Indeed, the 
condition number of B _1 A is order N 2 as opposed to the condition number of A which is order N 3 . 
This will be also the case for other problems and has already been noticed by [H. ] in a different 
context for an application to conjugate gradient algorithms. This can be viewed as an inverse 
inequality or a Jackson type one since the weighted L 2 -type norm associated to the bilinear form B 
is compared to the H^Aj-norm as follows 

(111.14) V <p € (P n (A) fl H^(A) , ||<p|| 2 A ^N(N-1)b(cp,tp) . 

Secondly, the other property that is very important in the multigrid algorithm is that the factor k 
as defined in Theorem 11.1 is also bounded by two constants independent of N. Indeed, we can state 
V ip € such that a^ GL (<p ,v) = 0 for any v in x" /2 
we have 

(111.15) b(ip ,ip) < k 1/2 a(ip ,(p) , 
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where the constant k is bounded by 4(N- 1 )/(N + 2) as can be derived from lemma 1 1 1.2 . 

Let us generalize now the results obtained we have obtained in this very simple situation to 
the case of a multielement discretization. 


i 


We begin here also bg analyzing the properties of the diagonal of the stiffness matrix A. We 
immediately note that there are two kinds of diagonal elements in this matrix : those that 
correspond to internal points, i.e. that involve the scalar product 
(IV. 1 ) A., = LL <?k,n [dhi/dx dh./dx]U nk ) 

with i □ 1 , and those that correspond to interface points, i.e,, that involve the scalar 

product (IV. 1 ) for i = 0 or N. The first kind of diagonal elements are the same as those involved in 
the previous section. Indeed, the corresponding Lagrangian interpolants vanish at the interfaces 
and also on any subinterval that does not contain the point ^ ; therefore, from Lemma 1 1 1. 1 , we 
can state that 

V i, i = 1 ,..,N-1 , V? ,i = 1 A u = (2/bp q. N(N+ 1 )/[3 ( 1 -(J;.) 2 )] , 
or again, thanks to a simple change of variable 


(IV. 2) 


,« = 1 >K , A ii< = o w N(N + 0/[3(^ r a < )(a t+1 -t; w )] 


For the interface terms we have 


Lemma IV. 1 : For i = 0 and any Z = 2....K we have 
(IV. 3) A., = [(b^ 1 ) -1 +(b < )- 1 ] (N 2 + N + l)/3 . 

P roof : As already used, the exactness of the Gauss- Lobatto formula gives, for i = 0 and Z = 2 ... ,K 
( or i = N and£ = 1 ,..,N- 1 , due to the direct stiffness summation) 

(IV. 4) A it = f A<-1 [dh N ,_/dx] 2 U) dt + | A< [dh 0 /dx] 2 (O dL . 

A simple change of variables and the use of the symmetry of h N and h 0 yields 
(IV. 5) \ t = [(2/b,_ 1 )+(2/b < )]f A [dh N /dx] 2 (t;) dt; . 

Let us compute the integral on the right-hand side of this equation. First we have 

(x+i) TT j N _ i (x-O 


h M (x) = 


zn,:, o-y 


F rom (111.13) we then get that 

h N (x)n (x+1)Ln(x)/ 2L'(1) , 
so that, after integration by parts, we obtain 
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} A [dh N /dx] 2 (x) dx = - (2 L^O))" 2 [j A (x + 1) L^(x) (d 2 /dx 2 )((x + 1) L^(x)) dx 

- 2 L^( 1 ) (d/dx)((x+ 1 ) L^(x))( 1 ) ] ; 

here again, the use of the exactness of the Gauss-Lobatto formula to compute the integral on the 


right-hand side yields 

(IV. 6) j A [dh N /dx] 2 (x) dx = - (2 L^( I))' 1 [(d 2 /dx 2 )((x+ 1 ) L^(x))( 1 ) p N 

- (d/dx)((x+ 1 ) L^(x))( 1 ) ] 


It is an easy matter to note that 


(d/dx) ( (x + 1 ) L’(x)) □ (x+ 1 ) L” (x) + L N (x) , 
(d 2 /dx 2 )((x+ 1 ) L^(x) ) o (x+ 1 ) Ln" (x) + 2 (x) , 

from (111.12) writen in the form 

( 1 -x 2 ) L„(x) - 2xL'(x) + N(N+1) L n (x) qO ; 
we derive that 


L^( 1 ) = N(N+ 1 )/2 , 

L^d) = (N- 1 )N(N+ 1 )(N + 2)/8 , 

L'" (1) = (N — 2)(N— 1 )N(N+ 1 )(N + 2)(N + 3)/48 ; 
this gives 

(d/dx)((x+ 1 ) L' (x))( 1 ) = N 2 (N+ 1 ) 2 /4 , 

(d 2 /dx 2 )((x+ 1 ) L^(x))( 1 ) = (N-1)N 2 (N+1) 2 (N + 2)/24 
Plugging this in (IV. 6) and using the relation (see [D.R.]) 

(IV. 7) q n = 2/(N+1)N , 
we derive 

(IV. 8) } A [dh N /dx] 2 (x) dx = (N 2 + N + 1 )/6. 

The lemma follows then from (IV. 5). 


N 

From (IV. 2) (IV. 3) and (IV. 7), we derive that the bilinear form b h GL that defines the 
smoother and is associated to the normalized diagonal part of A is proportional to the bilinear form 

N M 

b h GL defined for any <p and >j/ in X h as follows (remind (1 1 1.8)) 

(IV. 9) ^ >GL (<p,i|/) = b/ip.v) , 


where 
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6 ,(<p.«i») = (LflJ e w <ptt w ) 4 'U i _,) [U u - a P (a ^i-^P ] _1 

+ Po^W^o.P [ (n2+ n + 0/ (bp 2 ] + q m <pU m ) 4^,) [(N 2 + n + 1 )/ (bp 2 ] } . 

~N 

In the two next lemmas, we shall analyze the eigenvalue problem between a and b h CL that 

N ~N 

will allow first to estimate the normalization factor between b h GL and b h GL . 

Lemma IV. 2 : For any ip in xj , we have 

(IV. 10) j A tp , 2 (x) dx < (4N 3 + 2N 2 +3N-1)N/3(N 2 + N+1) 6 * 6l (tp ,tp) . 

Proof : To any (p in xj ) 1 let us associate the element <p 0 defined as follows 
(IV. 11) 4*o = 4* — ^,-2 4*( a ,) [^o ,, + h N ,,_i ] > 

so that ip 0 and 4 ) coincide on the interior collocation points. Using the inequality 
(a+b ) 2 < ( 1 +oc _1 ) a 2 + ( 1 +oc) b 2 , 
we deduce from (IV. 1 1 ) , 

(IV. 12 ) J A< 4 *' 2 (x) dx ^ ( 1 +a _1 ) J A| 4*o 2 (x) dx + ( 1 +a) J A< [4>(a,)ho ( + 4 *(a ?+ 1 )h^ ,] 2 (x) dx. 

It is an easy matter to note that the restriction of <p 0 to any A* belongs to IP N (Ap fl HpA,) so 
that the lemma 1 1 1.2 and a simple change of variables yields 

l A( 4*o 2 (x) dx< N(N- 1 ) Z-ll Pi., 4*o^i,,)[(C u - ap^+b,-^,)] -1 
< N(N- 1 ) L-l I 1 Pi,, 4* 2 (^,)[(^,- apca^+b,-^)]- 1 . 

Besides, from (IV. 8 ) , we derive that 

Ja, [h 0 ,,] 2 (x) dx = J A , [hp t ] 2 (x) dx = (N 2 + N + 1 )/3b, , 
as following the same lines as in the proof of lemma IV. 1 ; we get, for any N ^ 2 
1a, ho,, (x) dx ^ 1 /3bf . 

From (IV. 1 2) we then deduce 

f A ^' 2 (x) dx < (1 + a- 1 )N(N-1) L":\ Pj,, 4> 2 Ui,,)[ (K ir a ,)(a,+b ,-^,)]” 1 

+ ( 1 +a)[(N 2 +N+2)/3b ? ][ l 4 *(a ? )| 2 + l 4 *(a {+1 )| 2 ], 
choosing now oc = 3(N- 1 )(N 2 + N+ 1 )/(N + 1 )(N 2 +N+2); it follows from (IV. 7) and (IV. 9) that 
J A 4>' 2 (x) dx < (4N 3 + 2N 2 +3N-1)N/3(N 2 + N+1) b^ GL (4*.4*) , 
and the lemma follows. 

Remark IV. 1 : The estimate in the previous lemma provides a less precise characterization of the 
eigenvalue problem than the one we could get in lemma 1 1 1 . 1 , but we note that the highest 
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eigenvalue involves the same asymptotic order as in the previous section and this will be enough 
for our purpose. The important fact is that the result is independent of K and of the ratios between 
the various subinterval length b^ . As we shall see in what follows, this will result in a multigrid 
algorithm that will work as well in any case of number of subelements. Besides, note that the 
smallest eigenvalue of problem (11.8) scales like K" 2 (independently of N) such that the condition 
number of B _1 A behaves like (KN) 2 , in accordance with the finite element equivalent (when N is 
order 1 ) and proves that the conjugate gradient algorithm , when preconditioned by B , has a rate of 
convergence which behaves like 1 - c/KN and depends on both K and N ! ! ! This is of importance when 
we compare the preconditioned conjugate gradient with the multigrid algorithm. We refer to the 
thesis of E. M. RONQUIST [R] for numerical evidences. 


N 

It is an easy matter to derive from lemma IV. 2 that the normalized form b h GL is defined 

N o 

from b h 0L by multiplication by a factor of order ( 4/3) N . The other property that is important 
for the analysis of the multigrid algorithm deals with the space X^ /2i of those elements ip of X^ 1 
that verify 

(IV. 13) V t € X h N/2 ,a h N CL (4),ip)nO . 

Lemma IV. 3 : For any ip in X^ /2i we have 

(IV. 1 4) GL ((p ,(p) < 4(4N 3 + 2N 2 + 3N- 1 )/3(N 2 + N+ 1 )(N + 2) a(ip ,<p) . 


Proof : Let <p belong to X^ 21 , and let us define for £ = 2,..,K the element of X 

€ X^ , Vk = 1 K+ 1 , ^(a k ) = 8 <k . 

It is an easy matter to compute that 

1 /b ( _ } , for x in , 

<p^ (x) = - 1 /b, , for x in A ( , 

0 , for x in A k , k ^ £ and k ^ £-1 . 

Using this function in (IV. 1 3) , we derive that 

V 8 = 2 ,.. ,K , (1/b,_,)((p(a,) - <p(a, _,))+( 1/bp(<p(ap - <p(a t+1 )) = 0 ; 
recalling now that ip(a 0 ) □ ip(a K+1 ) = 0, we deduce that in fact <p vanishes at any 
ip ]A ^ is an element of [P N (A f ) D H G ( A f ). The use of (111.15) over each A ( proves 
b^ 1 CL ((P> l P) = b_ e ( ip ,cp ) ^ 4/N(N + 2) a(<p,<p) . 


N/2 

h 


by 


interface so that 
that 
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N 

and the lemma follows by recalling the normalization factor of b h GL . 

As a simple consequence of the previous lemmas and (11.1 0(11.14) we derive that the 
multigrid algorithm converges toward the numerical solution with a speed independent of both N 
and K, indeed, we have 

Theorem IV. 1 : The muJtigrid algorithm based on the Jacobi preconditioner converges, and at any 
V-cycle with m smoothing as detailed in section II , the rate of convergence is given by 
J A e; 2 (x) dx < * j A e G 2 (x) dx , 

where 

t = [(k- 1 )/«] 2m , for m < (k-1 ) , 
and 

X □ (Kf(m , 1 )) 2 , for m>(K-1), 

with 

K □ 4(4N 3 + 2N 2 + 3N-1)/3(N 2 + N + 1)(N + 2) = (16/3) (1 +0(N -1 )) . 

Remark IV. 2 : Let us note that the convergence rate that we have theoretically obtained is 
independent of K and the sizes of the subintervals and does not deteriorate when N increases; this is 
in accordance with what is numerically observed in part one of this paper [R.P.], However, the 
rather rough estimate we gave for the highest eigenvalue in (IV. 10) provides a rather too high 
estimate for the convergence parameter ( close to 0.8 1 when N is large enough ). Note also that the 
optimal choice of parameter 4 < m < 5 is close to what can be observed numerically. As in the 
analysis of the multigrid algorithm when a finite element method is used, this optimal parameter 
is lower than the actual one. By using negative order of the norms defined from b and a, we could fit 
more closely to the experiments for this last result. This would be obtained, however, at the price 
of a much more technical proof and would not really be worthy since the main conclusion is the 
independence of the convergence rate with respect to the parameters of the discretization, 
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V. Anal us i s of the convergence of the multiarid algorithm in the 
case of non constant coefficients. 

The previous chapter was devoted to the analysis of the multigrid algorithm when applied to 
the very simple equation -u xx = f. Of interest of course is to know that the same conclusions hold 
true in more complicated situations. We shall extend here the analysis to the case of the equation 
(V. 1 ) -(cxu x ) x = f , 

where a is a function of x such that there exists two constants ex' and a + with 
(V.2) Vxe A , 0 < ex' < cx(x) < a + , 

and also such that a is in the Sobolev space W 2 - 00 (A). We shall assume here that the domain is not 
decomposed into subdomains, leaving this analysis to a forthcoming paper as well as the analysis of 
the multidimensional case. The analysis provided here is inspired from the reference [B.V.], 

It is standard to note here that problem (V. 1 ) can also be stated in a variational formulation 
like (11.19) with now a defined as follows 
(V.3) V ip , 4 / 6 Hq( A) , a(<p ,40 = J A cx(x) ip’(x) 4 /'(x) dx . 

It is rather well known also that the general spectral(element) discretization of the equation 
consists in first choosing a discretization parameter n € IN*, then : Find u n € X n = (P n (A) f| Hq(A) 
such that 

(V.4) V v n € X n , a n (u n ,v n ) = (f,v n ) n , 
where the discrete form a n is defined here by 
(V.5) V ip , 41 € X n , a n (ip , 4 /) = (oc cp x , ip x ) n . 

The nested spaces for the multigrid algorithm are exactly the same as in the previous sections, i.e. , 

6 At 1 = X N/2 and cM 2 = X N , and the strategy also based on the use of the Jacobi preconditioner as a 

smoother. The corresponding bilinear form b N is deduced from the following 6 N after normalization 

of the maximum eigenvalue of problem (11.7): 6 N is defined by 

(Y. 6 ) V 4 ) , 4 < € X N , 6 N ((p,i|/) = L „ =0 Ann<P (£„)¥(£„) Pn - 

where A nn corresponds to the diagonal entry of the stiffness matrix A, equal to 

A . J^kN k Nn 

and we recall that h„ is the Lagrangian interpolant at the point ^ .From hypothesis (V.2) and the 
exactness of the Gauss-Lobatto quadrature formula, it is simple to derive that 
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(V.7) V cp € X N , a N ((p,ip)< a + a(tp ,tp) < a + N(N- 1 ) 6 (41 ,<p) 

< ( 3 a + /oc”)(N - 1 )/(N+ 1 ) 6 N (<p,<p), 

where we recall that B( tp ,»|i) = J A ip(x) \p(x)/( 1 -x 2 ) dx. As a result, we derive that the 

normalization of b N satifies 

(V. 8) b N = (3oc + /oc‘)(N - 1 )/(N + 1 ) 6 N . 

From the eigenvalue problem : Find ^ in X N and X in [R + such that 
(V. 9 ) V v € X N , a N (^,v) = X b N (^,v) 

that possess N- 1 eigenvalues 0 < X 1 X 2 < X N-1 < 1 , the eigenvectors of which are chosen 
normalized with respect to the norm derived from b N we define, as in (II. 1 0) the ||| . ||| 0N norms for 
any 0 € [R by 

(v. 1 0 ) V V € X" , III V ||| e 2 N = L":! xj b"(y,f„) 2 . 

Here the bilinear ^orms associated to the definition of the system and to the smoothing depend on N. 

The definitions of the smoothing operator^ and the correction operator G N have to be precised. Let 

us do this in the simple case where only two grids are used. Instead of (11.3), the smoothing 

procedure consists in : Find in X N such that 

(V. I 1 ) V v € X N , b N (A-«P,v) = g N (v)-a N (ip,v) , 

while the correction procedure consists in : Find in X N/2 such that 

(V. 1 2) V v € X N/2 , a N/2 (ip,v) □ g N (v)-a N ((p ,v) , 

and define 0% = ip + ^p. Then as explained in the general case, during a V-cycle with m/2 
smoothings down and up the error e° is changed in e 1 as follows 
(V. 1 3) e 1 = (£> N ) m/2 G n ($ N ) m/2 e° . 

Before entering in the details of the analysis of the decay rate of the error, let us state some 
results of general interest. First of, let x be in IP°(A) be given, r(x) and r N (x) in 1P^(A) be such 
that 

(V. 1 4) V 4 / € P n °(A) , b(r(x).»l») = b N (r N (xU) = a N (x,4-) , 

where b is defined in (1 1 1. 1 3). It is a simple consequence of (V.9) to derive that 

III r N (x) ||| 2 n = b N (r N (x),r N (x)) = a N ( y ,r N ( x)) < ||| x lll 2 , N III r N (x) III 0 .n > 
nfj( X) |||q jg ^ III X III2 ,N ' 


whence 


- 22 - 


Besides, from (Y.7) and (V. 1 4) , we have 

b(r( X ),r( X )) - b N (r N ( X ) ,r(x)) [b N (r N ( X ) ,r N ( X )) b N (r( X ) ,r( X ))] 1/2 

^ (a +2 /oO ,/2 [ b N ( r N ( X ) ,r N ( X )) b(r( X ) ,r( X ))] 1/2 ; 
finally, we derive that the solution r( X ) of (V. 1 4) has the following stability property 
(V.15) b(r( X ),r( X )) < c||| X ||| 2N . 

Then let us state some approximation result the proof of which will be presented in the appendix. 
Lemma V. 1 : Let p defined over A and such that 

(V. 16) | A p 2 (x)/( 1 -x 2 ) dx < oo , 

and n be a positive integer. The solutions <p(p) € Hq(A), tp n (p) € P°(A) defined by 
(V. 1 7) V v € Hq(A) , a(<p(p),v) = b(p,v) , 

(V. 1 8) V v € (P°(A) , a n (ip n (p),v) = b(p,v) ; 

then the following approximation results holds 
(V. 1 9) I <p(p) - <p n (p) I? < c n' 2 B(p,p) , 

(V.20) 6(<p(p) - 4» n (p) ,tp(p) - <p n (p)) < cn -4 b(p,p) . 

Let us denote by e the term ($ N ) m/2 e 0 ; derive now as in (1 1. 1 2) that, 

(V.2 1 ) III e 1 III, N < f(m/2,1 /2) ||| 6 N e ||| 0N , 

( V . 2 2 ) |||e||| 2N <ftm/2,l/2)|||e 0 ||| 1N . 

From the definition of (3 N , we derive that, for any X in IP^(A) 

III e N X ||| 0 N = III X - x IIIo.n < HI X - N(N- 1 )cp(r( X )) ||| 0 N + III X - N(N— 1 )ip(r( X )) ||| 0 N ; 

hence 

III e N X ||| 0N < (a +2 /a')[ b( X - N(N-1 )<p(t(x), X - N(N- 1 )tp(r( X )) 1/2 

+ b(x - N(N- 1 )cp(r( x) , X - N(N-1 )<p(r( X )) 1/2 ] . 
With the previous notations (see (V. 1 4) and (V. 1 8)) and recalling that the normalization factor 
between the forms b and b is N(N-l) (see (111.14)), it is simple to derive that 
X = N(N- 1 )<p N (r( X )), while "x (defined in (V. 1 2)) satisfies^ = N(N- 1 )<p N/2 (r( X )); we derive 
III e" X ||| 0N < (oc +2 /oc”)(N(N- 1 ))[ b(ip N (r( X )) - <P (r( X ), «p N (r( X )) - <p(r( X )) 1/2 

+ b( <P N/2 (r( X )) - <p(r( X ), tp N/2 (r( X )) - ip(r( X )) 1/2 ] , 

applying next (V.20) for n = N and n = N/2, then (V. 1 5) , we obtain 
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Lemma V.2 : There exists a constant c independent of both N and K such that 
(V.23) V x 6 Pj(A) , ie N li ||| 0H ^c]||xlll 2 .H ■ 

Plugging (V,23) in (V.2 1 ) and using (V.22) give now the following 

Theorem V. 1 : The multigrid algorithm based on the Jacobi preconditioner converges, and at any 
V-cycle with m smoothings as detailed in section II , the rate of convergence is given by 
(V.24) } A e; 2 (x) dx ^ c /tm,l ) J A eQ 2 (x) dx . 

Remark V. 1 : The decay rate of multigrid algorithm is independent of N, and it is interesting to 
note that its asymptotic behavior is 0( 1 /m). 

Remark V.2 : The analysis of the case of non constant coefficients requires some regularity (i.e. , 
W 2,0 ° ) of the coefficient ( this is required for lemma V. 1 ). We do not know to what extend the lack 
of regularity of oc deteriorates the actual convergence. Note however that the local regularity is 
just required, i.e., oc can be non smooth through the interfaces. The multigrid procedure proposed 
in [ZWH 1.2], and that is based on another approach of smoothing, is numerically proved to be very 
robust with respect to the irregularity of oc [*]. 

A. APPENDIX 

The main purpose of this appendix is to provide the proof of lemma V. 1 . Let us first recall 
the following result of the approximation theory (see [D]) 

Lemma A. 1 : Let n M denote the orthogonal projection operator from L 2 (A) onto IP M (A). The 
following approximation results hold 

(A. 1 ) V u € H m (A) , II u - n M u || 0 < c(m) M _m [ J A u (m)2 (x) ( 1 -x 2 ) m dx] 1/2 . 

Proof : Let us recall that the Legendre polynomials constitute a total system of orthogonal 

functions of L 2 (A) that verifies 

(A-2) J A L„(x)L m (x)dx = 8„,2/(2n t l) . 


* T.A.ZANG- personnal communication. 
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Let us use this basis to span the function u; we arrive at 
^ =0 ’ 

so that the L 2 (A)-projection of u onto IP M (A) is equal to 

*•* - ^n :0 ’ 

and the error 

(A.3) u-n M u=L n ~ M+1 a n L n . 

Let us recall now that the Legendre polynomials satisfy the following relation 
( ( 1 -x 2 ) L‘ ) = - n(n+ 1 ) L n . 

From (A. 2) we conclude that 

(A. 4) J A L n (x) L^x) ( 1 -x 2 ) dx = 8 nm 2n(n+1 )/(2n+ 1 ) . 

From (A. 2) and (A.3) we get 

l|u-n M u|| 0 = L n ~ M+1 2 a 2 /(2n+ 1 ) , 
while, from (A. 4) , we derive that 

j A u l2 (x) ( 1 -x 2 ) dx □ L~o 2n(n+ 1 ) a 2 /(2n+ 1 ) , 
and (A. 1 ) is then just a simple consequence of these two equalities in the case m = 1 . The general 
case is handled by recursion. 

As a consequence, we derive that 

Corollaru A. 1 : Let ti^ denote the orthogonal projection operator from Hq(A) onto tPjjJ(A). 
The following approximation results hold 

(A. 5) V u € H m (A) fl Hq ( A ) , || u - t^u ||, < c(m) M _m [ } A u (m+1)2 (x) ( 1 -x 2 ) m dx] 1/2 . 
Proof : It has already been noted that for any u in Hq(A) 

[tT|J|U](x) = \*_ y n M _ 1 (u‘)(t) dt , 
so that (A. 5) is an easy consequence of lemma A. 1 . 

Proof of lemma V.1 : It is an easy matter to note that the solution <p(p) of problem (V.17) 
satisfies 

V x € A , - [(oc[<p(q)] x )J(x) = p(x)/( 1 -x 2 ) . 

From (V.2) and (V. 1 6) , we derive 

(A. 6) J A [<p(p)] 2 (x) dx + J A [<p(g)] 2 x (x) (1 -x 2 ) dx < c } A p 2 (x)/( 1-x 2 ) dx < c B(q,q) . 

From corollary A. 1 , we know that there exists an element i|/ n of P^CA), such that 
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I <p(p) - >V n Ii < c n' 1 J A [tp(Q)]J x (x) ( 1 -x 2 ) dx , 
and 

IvJ, <1 <p(e) I, , 

so that 

(A. 7) | «p(p) - M# n h < c n" 1 b(e , q) . 

From the ell i pticity of the form a 11 , we derive that there exists a constant §, such that 

& I - Y n I? < a n («p n (p) - «|» n , «p n (o) - V n ) i 

using (V. 1 7) and (V. 1 8) we derive that 

£ I <P n (p) - I? < a(<p(p), «P„<e) - ¥ n ) ~ a n (^ n . <P n (p) - Y n ) . 

or again 

(A. 8) & | <p n (e) - v|/ n | 2 < a(<p(p) - n/ n , tp n ( q) - v n ) + a(i|/ n , «p n (o) - <|/ n ) - a n (i|/ n , <p n (p) - y n ) • 

Let us note now that from the exactness of the Gauss Lobatto quadrature formula we get that 
(A.9) v OC n € P n/2 (A), (*y n ,(ip n (p) - <|» n )') = .<«P n <Q) - Yn)') n> GL > 

so that 

a(^ n , <P n (p)-Y n ) - a n (^ n , <P n (e)-Y n ) = ((a-oc n )*i ,(<P n (p)-»l' n )) - ((<x-oc n )^ ,( <P n ( Q)-^ n ) ) n • 
and using now (A. 8) yields 

(A. 10) |(p n (p) - <|i„ |, < c( |<p(p) - I, + I OC - « n l|_~(A) I h ) • 

It is standard to note that there exists an element a n such that 
l 0l - a nlL oo (A) <cn " k ^j k =ol« (j) lL“(A) ; 
hence, from (A. 1 0) and (A. 7) , we derive that 

I <p(p) - <P n (Q) I, < c n' 1 (b(p ,p) 1/2 + (Lj l0 1 (X<J> 1 l°°( A)^ I ^ li ^ - 

and (V. 1 9) is an easy consequence of (A.6). 

Let us turn now to the proof of (V.20). We shall use a standard duality technique and first 
define an element x in Hq(A) as follows 
(A. 11) Vv€Hj(A) , a(x,v) = b(<p(p) - <p n (p), v) , 

1 1 is an easy matter to derive 

b(tp(p) - <p n (p),<p(p) - <p n (p)) = a(x, <p(p) - <p n (p)) 

= a(x - TT^ x , <p(p) - <p„(p)) + a(n,J /2 x , <p(p) - «p„(p)) 

= - nl/2 X - <p(p) - + (a " - a)(n rl/2 X , <P„(p)) • 

Let us examine the last term on the right-hand side of this last equality; using (A.9) one more time 
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we write 

(a n - a)(n r J /2 x , <P n (p)) = ((oc - a n )( n' n/2 x)', <P n (p)') n 

+ ((oc - oc n )( n^ /2 x)', tp n (Q)‘) . 

so that 

b(«p(p) - <p n (p),<p(p) - <P n (p)) = a(x - n„/ 2 x.<p(p) - <P n (p)) + ((oc - oc n )(n^ /2 x)', <p n (p)’) n 

+ (Coe - oc n )( n^ /2 x)', tp n (Q)’) ; 

hence 

b(«p(p) - <p n (p),<p(p) - <p n (p)) <c| X - K /2 X li I <P(Q) - ¥ n (p) It 

+ c n' 2 Lj 2 =0 1 oc (j) l L «(A) I n n/2 X It I <P n (p) li > 
< c| X - Tt,|/ 2 X I, I <p(p) - <p n (p) I, 

+ C n’ 2 Lj 2 :0 I 0C (j) l L °° (A) I n n/ 2 X I, U n (p)li . 
Using now corollary A. 1 and (A. 1 1 ), we derive (after bounding I oc (j) by a constant) 

[b(tp(p) - ip n (p),ip(e) - ip n (p))] 1/2 < c (n _1 |ip(p) -ip n (p)|, + n' 2 | ip n (p) |, ) ; 
thanks to the stability of ip n (p) with respect to p, we derive (V.20) from (V. 1 9), 
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